*
* Table T01 T12 Double Robust
*
set varabbrev off, permanently
clear all
use epfipc2006-2019plus.dta
*
* Total expenditure
*
global listtotal gasto1-gasto12 year2-year13 ccaa2-ccaa18 income squaredm othehouse ///
       ttenure2-ttenure6 tsize2-tsize5 pdensity2-pdensity3 residshort nmembers nfemales ///
	   nforeign ndworkers n5_15 n16_24 n25_34 n35_64 n65_79 n80omas nstudies2-nstudies4 ///
	   nactive nincome labinc peninc uneinc capinc agemb femalemb foreignmb unemployedmb inactivemb
*
*
foreach var1 in T01 T12 {
	preserve
	
	dis in blue "`var1'"
    keep if post==0 & `var1'!=.
	* propensity score estimation with survey weights
    quietly logit `var1' $listtotal  [pweight=factor]
	predict double ps if e(sample)==1, pr
	drop if ps==.
    gen w=(`var1'+(((1-`var1')*ps)/(1-ps)))
	* Weights used for ATET estimation accounting for survey weights
	gen estweight=w*factor
	* imposing common support
	quietly sum ps if `var1'==1
    scalar pmintreat=r(min)
	drop if ps<pmintreat & `var1'==0
	quietly sum ps if `var1'==0
    scalar pmaxuntreat=r(max)
	drop if ps>pmaxuntreat & `var1'==1
	* estimating ATETs
	foreach v of varlist gasto1-gasto12 year2-year13 ccaa2-ccaa18 income squaredm othehouse ///
       ttenure2-ttenure6 tsize2-tsize5 pdensity2-pdensity3 residshort nmembers nfemales ///
	   nforeign ndworkers n5_15 n16_24 n25_34 n35_64 n65_79 n80omas nstudies2-nstudies4 ///
	   nactive nincome labinc peninc uneinc capinc agemb femalemb foreignmb unemployedmb inactivemb	{
	   quietly sum `v'
	   generate dm_`v' = `v' - r(mean)
	   gen inter`v' = dm_`v' * `var1' 
	}
	*
	global demeaned dm_gasto1-dm_gasto12 dm_year2-dm_year13 dm_ccaa2-dm_ccaa18 dm_income dm_squaredm ///
	dm_othehouse dm_ttenure2-dm_ttenure6 dm_tsize2-dm_tsize5 dm_pdensity2-dm_pdensity3 dm_residshort ///
	dm_nmembers dm_nfemales dm_nforeign dm_ndworkers dm_n5_15 dm_n16_24 dm_n25_34 dm_n35_64 dm_n65_79 dm_n80omas ///
	dm_nstudies2-dm_nstudies4 dm_nactive dm_nincome dm_labinc dm_peninc dm_uneinc dm_capinc dm_agemb dm_femalemb ///
	dm_foreignmb dm_unemployedmb dm_inactivemb 
	*
	global interaction intergasto1-intergasto12 interyear2-interyear13 interccaa2-interccaa18 interincome /// 
	intersquaredm interothehouse interttenure2-interttenure6 intertsize2-intertsize5 interpdensity2-interpdensity3 ///
	interresidshort internmembers internfemales internforeign interndworkers intern5_15 intern16_24 intern25_34 ///
	intern35_64 intern65_79 intern80omas internstudies2-internstudies4 internactive internincome interlabinc ///
	interpeninc interuneinc intercapinc interagemb interfemalemb interforeignmb interunemployedmb interinactivemb 
	*
	* needed to order regressors
	order dm_* inter*, sequential
	*
	forvalues cat=1/12 {
	    display in red "`: var la gasto`cat''"
		* estimation using survey weights
	    reg dgasto`cat' `var1' $demeaned $interaction [aweight=estweight], r
		if `cat'==1 {
			count if `var1'==1 & e(sample)==1
			local ntreated=r(N) 
			count if `var1'==0 & e(sample)==1
			local nuntreated=r(N)
			outreg2 using table`var1', excel dec(2)  adds(number_treated, `ntreated', ///
			number_untreated, `nuntreated') replace
		}
		else {
			* count number of treated in ATET estimation
			count if `var1'==1 & e(sample)==1
			local ntreated=r(N) 
			count if `var1'==0 & e(sample)==1
			local nuntreated=r(N)
			outreg2 using table`var1', excel dec(2)  adds(number_treated, `ntreated', ///
			number_untreated, `nuntreated')
		}
	}
	* auxiliary regressions
	quietly reg dgastot $demeaned $interaction [aweight=estweight]
	predict udg, residuals
	quietly reg dincome $demeaned $interaction [aweight=estweight]
	predict udy, residuals
	quietly reg `var1' $demeaned $interaction [aweight=estweight]
	predict ut, residuals
	
	gmm (eq1: udg - {beta1}*ut) /// 
    (eq2: udy - {beta2}*ut) [pweight=estweight], ///
	instruments(eq1: ut) instruments(eq2: ut) winitial(identity) twostep
	
	nlcom (/beta1)/(/beta2)
	nlcom (/beta1)/(/beta2)-1
	display "Number of treated"
	count if `var1'==1 & e(sample)==1
	display "Number of untreated"
	count if `var1'==0 & e(sample)==1
	
	restore
}
*
* Expenditure per equivalent adult
*
global listtotalaeq gastoaeq1-gastoaeq12 year2-year13 ccaa2-ccaa18 income squaredm othehouse ///
       ttenure2-ttenure6 tsize2-tsize5 pdensity2-pdensity3 residshort nmembers nfemales ///
	   nforeign ndworkers n5_15 n16_24 n25_34 n35_64 n65_79 n80omas nstudies2-nstudies4 ///
	   nactive nincome labinc peninc uneinc capinc agemb femalemb foreignmb unemployedmb inactivemb
*
foreach var1 in T01 T12 {
	preserve
	dis in blue "`var1'"
    keep if post==0 & `var1'!=.
	* propensity score estimation with survey weights
    quietly logit `var1' $listtotalaeq  [pweight=factor]
	predict double ps if e(sample)==1, pr
	drop if ps==.
    gen w=(`var1'+(((1-`var1')*ps)/(1-ps)))
	* Weights used for ATET estimation accounting for survey weights
	gen estweight=w*factor
	* imposing common support
	quietly sum ps if `var1'==1
    scalar pmintreat=r(min)
	drop if ps<pmintreat & `var1'==0
	quietly sum ps if `var1'==0
    scalar pmaxuntreat=r(max)
	drop if ps>pmaxuntreat & `var1'==1
	* estimating ATETs
	foreach v of varlist gastoaeq1-gastoaeq12 year2-year13 ccaa2-ccaa18 income squaredm othehouse ///
       ttenure2-ttenure6 tsize2-tsize5 pdensity2-pdensity3 residshort nmembers nfemales ///
	   nforeign ndworkers n5_15 n16_24 n25_34 n35_64 n65_79 n80omas nstudies2-nstudies4 ///
	   nactive nincome labinc peninc uneinc capinc agemb femalemb foreignmb unemployedmb inactivemb	{
	   quietly sum `v'
	   generate dm_`v' = `v' - r(mean)
	   gen inter`v' = dm_`v' * `var1' 
	}
	* 
	global demeaned dm_gastoaeq1-dm_gastoaeq12 dm_year2-dm_year13 dm_ccaa2-dm_ccaa18 dm_income dm_squaredm ///
	dm_othehouse dm_ttenure2-dm_ttenure6 dm_tsize2-dm_tsize5 dm_pdensity2-dm_pdensity3 dm_residshort ///
	dm_nmembers dm_nfemales dm_nforeign dm_ndworkers dm_n5_15 dm_n16_24 dm_n25_34 dm_n35_64 dm_n65_79 dm_n80omas ///
	dm_nstudies2-dm_nstudies4 dm_nactive dm_nincome dm_labinc dm_peninc dm_uneinc dm_capinc dm_agemb dm_femalemb ///
	dm_foreignmb dm_unemployedmb dm_inactivemb 
	*
	global interaction intergastoaeq1-intergastoaeq12 interyear2-interyear13 interccaa2-interccaa18 interincome /// 
	intersquaredm interothehouse interttenure2-interttenure6 intertsize2-intertsize5 interpdensity2-interpdensity3 ///
	interresidshort internmembers internfemales internforeign interndworkers intern5_15 intern16_24 intern25_34 ///
	intern35_64 intern65_79 intern80omas internstudies2-internstudies4 internactive internincome interlabinc ///
	interpeninc interuneinc intercapinc interagemb interfemalemb interforeignmb interunemployedmb interinactivemb 
	*
	order dm_* inter*, sequential
	forvalues cat=1/12 {
	    display in red "`: var la gasto`cat''"
		* estimation using survey weights
	    reg dgastoaeq`cat' `var1' $demeaned $interaction [aweight=estweight], r
		if `cat'==1 {
			count if `var1'==1 & e(sample)==1
			local ntreated=r(N) 
			count if `var1'==0 & e(sample)==1
			local nuntreated=r(N)
			*outreg2 using tableaeq`var1', excel dec(2)  adds(number_treated, `ntreated', ///
			number_untreated, `nuntreated') replace
		}
		else {
			* count number of treated in ATET estimation
			count if `var1'==1 & e(sample)==1
			local ntreated=r(N) 
			count if `var1'==0 & e(sample)==1
			local nuntreated=r(N)
			*outreg2 using tableaeq`var1', excel dec(2)  adds(number_treated, `ntreated', ///
			number_untreated, `nuntreated')
		}
	}	
	* auxiliary regressions
	quietly reg dgastotaeq $demeaned $interaction [aweight=estweight]
	predict udg, residuals
	quietly reg dincomeaeq $demeaned $interaction [aweight=estweight]
	predict udy, residuals
	quietly reg `var1' $demeaned $interaction [aweight=estweight]
	predict ut, residuals
	
	gmm (eq1: udg - {beta1}*ut) /// 
    (eq2: udy - {beta2}*ut) [pweight=estweight], ///
	instruments(eq1: ut) instruments(eq2: ut) winitial(identity) twostep
	
	nlcom (/beta1)/(/beta2)
	nlcom (/beta1)/(/beta2)-1
	display "Number of treated"
	count if `var1'==1 & e(sample)==1
	display "Number of untreated"
	count if `var1'==0 & e(sample)==1
	
	restore
}
*
* Budget shares
*
global listtotalbs share1-share12 year2-year13 ccaa2-ccaa18 income squaredm othehouse ///
       ttenure2-ttenure6 tsize2-tsize5 pdensity2-pdensity3 residshort nmembers nfemales ///
	   nforeign ndworkers n5_15 n16_24 n25_34 n35_64 n65_79 n80omas nstudies2-nstudies4 ///
	   nactive nincome labinc peninc uneinc capinc agemb femalemb foreignmb unemployedmb inactivemb
*
foreach var1 in T01 T12 {
	preserve
	dis in blue "`var1'"
    keep if post==0 & `var1'!=.
	* propensity score estimation with survey weights
    quietly logit `var1' $listtotalbs  [pweight=factor]
	predict double ps if e(sample)==1, pr
	drop if ps==.
    gen w=(`var1'+(((1-`var1')*ps)/(1-ps)))
	* Weights used for ATET estimation accounting for survey weights
	gen estweight=w*factor
	* imposing common support
	quietly sum ps if `var1'==1
    scalar pmintreat=r(min)
	drop if ps<pmintreat & `var1'==0
	quietly sum ps if `var1'==0
    scalar pmaxuntreat=r(max)
	drop if ps>pmaxuntreat & `var1'==1
	* estimating ATETs
	foreach v of varlist share1-share12 year2-year13 ccaa2-ccaa18 income squaredm othehouse ///
       ttenure2-ttenure6 tsize2-tsize5 pdensity2-pdensity3 residshort nmembers nfemales ///
	   nforeign ndworkers n5_15 n16_24 n25_34 n35_64 n65_79 n80omas nstudies2-nstudies4 ///
	   nactive nincome labinc peninc uneinc capinc agemb femalemb foreignmb unemployedmb inactivemb	{
	   quietly sum `v'
	   generate dm_`v' = `v' - r(mean)
	   gen inter`v' = dm_`v' * `var1' 
	}
	* 
	global demeaned dm_share1-dm_share12 dm_year2-dm_year13 dm_ccaa2-dm_ccaa18 dm_income dm_squaredm ///
	dm_othehouse dm_ttenure2-dm_ttenure6 dm_tsize2-dm_tsize5 dm_pdensity2-dm_pdensity3 dm_residshort ///
	dm_nmembers dm_nfemales dm_nforeign dm_ndworkers dm_n5_15 dm_n16_24 dm_n25_34 dm_n35_64 dm_n65_79 dm_n80omas ///
	dm_nstudies2-dm_nstudies4 dm_nactive dm_nincome dm_labinc dm_peninc dm_uneinc dm_capinc dm_agemb dm_femalemb ///
	dm_foreignmb dm_unemployedmb dm_inactivemb 
	*
	global interaction intershare1-intershare12 interyear2-interyear13 interccaa2-interccaa18 interincome /// 
	intersquaredm interothehouse interttenure2-interttenure6 intertsize2-intertsize5 interpdensity2-interpdensity3 ///
	interresidshort internmembers internfemales internforeign interndworkers intern5_15 intern16_24 intern25_34 ///
	intern35_64 intern65_79 intern80omas internstudies2-internstudies4 internactive internincome interlabinc ///
	interpeninc interuneinc intercapinc interagemb interfemalemb interforeignmb interunemployedmb interinactivemb 
	*
	order dm_* inter*, sequential
	forvalues cat=1/12 {
	    display in red "`: var la gasto`cat''"
		* estimation using survey weights
	    reg dshare`cat' `var1' $demeaned $interaction [aweight=estweight], r
		if `cat'==1 {
			count if `var1'==1 & e(sample)==1
			local ntreated=r(N) 
			count if `var1'==0 & e(sample)==1
			local nuntreated=r(N)
			*outreg2 using tablebs`var1', excel dec(4)  adds(number_treated, `ntreated', ///
			number_untreated, `nuntreated') replace
		}
		else {
			* count number of treated in ATET estimation
			count if `var1'==1 & e(sample)==1
			local ntreated=r(N) 
			count if `var1'==0 & e(sample)==1
			local nuntreated=r(N)
			*outreg2 using tablebs`var1', excel dec(4)  adds(number_treated, `ntreated', ///
			number_untreated, `nuntreated')
		}
	}	
    restore
}
/*
*
* Invoke OS command to merge table files
*
! cat tableT01.txt tableT12.txt tableaeqT01.txt tableaeqT12.txt tablebsT01.txt tablebsT12.txt > table_T01T12DR.csv
* Delete individual table files
! rm table*.txt	
! rm table*.xml
*/
* The End
